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ABSTRACT 

We compare the accuracy of various methods for determining the transfer of the 
diffuse Lyman continuum in H ii regions, by comparing them with a high-resolution 
discrete-ordinate integration. We use these results to suggest how, in multidimensional 
dynamical simulations, the diffuse field may be treated with acceptable accuracy with- 
out requiring detailed transport solutions. 

The angular distribution of the diffuse field derived from the numerical integra- 
tion provides insight into the likely effects of the diffuse field for various material 
distributions. 
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1 INTRODUCTION 

In this paper, we model the diffuse field structure of as- 
tropliysical Hll regions. Recombinations within these re- 
gions emit radiation in the various hydrogen line spectra and 
Balmer and higher continua which is an observed character- 
istic of them. Emission in the Lyman continuum, however, 
is energetic enough to ionize other hydrogen atoms, and so 
is trapped within the nebula. This radiation field is believed 
to have a significant fraction of the energy density of the 
direct ionizing continuum in some parts of the Hll region, 
and so it is important to model it correctly. This is particu- 
larly relevant when modelling the complex dynamics of the 
nebulae, or the emission from the internal features such as 
the tails of cometary globules in the Helix planetary nebula 
( O'Dell et al.|2007||. 



Ritzerveld ( 2005 1 has suggested that diffuse fields may 
be particularly important at the edges of H ll regions, in 
some regimes. This would seem to suggest that the diffuse 
field may have a more significant impact on their overall evo- 
lution than previously assumed. However, Ritzerveld uses a 
simple outward-only treatment of the diffuse radiation trans- 
port, and also assumes that the absorption coefficient is com- 
parable for diffuse and direct radiation fields, although it is 
acknowledged that in reality the direct photons are likely to 
have a harder spectrum and will thus be significantly more 
penetrating. 

To investigate the validity of these conclusions, in this 
paper we apply detailed discrete-ordinate angular integra- 
tion to investigate the validity of several approximate nu- 



* (c) Crown Copyright 2009/MoD 

t E-mail: robin.williams@awe.co.uk 

w.henney@astrosmo.unam.mx (WJH) 



(RJRW); 



merical transport schemes. To simplify the problem, we as- 
sume a pure-hydrogen nebula, without dust, and use a sim- 
ple two-frequency approximation to the radiation flow. 

With our more detailed modelling, we find that the dif- 
fuse field can indeed dominate for the situations Ritzerveld 
describes. However, where the diffuse field dominates, it will 
typically also be outwardly-beamed and therefore be indis- 
tinguishable from the direct field for most purposes. For 
most astrophysically relevant conditions, the usual on-the- 
spot approximation is shown to give accurate results through 
most of a spherical nebula, except for a region close to the 
star where it overestimates the diffuse field. 

We also assess the accuracy of the Eddington diffusion 
approximation for the diffuse field transfer, which may be a 
useful means of modelling diffuse transport effects in multi- 
dimensional simulations. 

In the present paper, we neglect the effects of dust and 
heavy elements in the Hll region. This is a reasonable as- 
sumption for the case of cosmological H ll regions; however, 
there is observational evidence for the importance of dust 
absorption within the ionized gas in galactic H ll regions 
( [Cesarsky et al.|2000|[Robberto etaLpO OS ). We briefly dis- 
cuss the impact dust might have on our results. 



2 GOVERNING EQUATIONS AND PHYSICAL 
CONDITIONS 

We use a simple model of H ll regions to study the diffuse 
field structure. We consider steady, spherically symmetric 
pure- hydrogen regions, with the radiation field in two fre- 
quency components: higher frequency radiation propagating 
directly from a central point source, and lower frequency dif- 
fuse radiation. We assume that the material within the H ii 
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region is maintained at a constant temperature in ionization 
equilibrium. 

In the absence of dust, the governing equations are then 



dLdir/dr 



-n(l - a:)ao-Ldir 



didii/ds = —n(l — x)aildi{ + -p-T^^x^ 



(1) 
(2) 



2 2 



n{l — x) 



aoLdiv 



+ ai / Idii{r,n')d^' (3) 



corresponding to the transfer of direct radiation from a point 
source, Ldir, the diffuse radiation field /dif (r, 9), and balance 
of ionization and recombination for ionization fraction x. 
Here L is the specific luminosity and / the radiation intensity 
specified in terms of the number of photons, rather than 
the radiation energy, as this is more natural in the present 
application. 

Integrating equations ([l|-([3| over volume, and angle 
where appropriate, and applying the boundary conditions 
LdiT = at r = 0, Jdif = as r — > oo, we find that 



Li, = I aBU^x^dV, 



(4) 



where cte ~ Qa — "i- This has a simple interpretation, in 
that the number of ionizing photons into the nebula balances 
the number of recombinations which destroy such photons 
(rather than re-emitting a Lyman continuum photon). 

If the ionization front is thin, as is typically the case, 
then a characteristic Stromgren radius for a spherical H ll 
region, rs, is given by 



(5) 



independent of diffuse radiation transport effects. 



From Osterbrock and Ferland (20061, we take values 
of the recombination coefficients at 10 K of oa = 4.18 x 
10"^^cm^s"^ and Qi = 1.58 x 10""cm^s"^ We assume 
that the ionization cross section for the diffuse radiation is 
the Lyman limit cross section ai — 6 x 10^^^ cm^, while 
for the direct photons, we take either ao = 1 x 10~^*cm^ 
(which allows for the largest likely effect of the higher aver- 
age frequency of the direct photons) or ao = 6 x 10~^* cm'^ 
(corresponding to the case which was primarily considered 
by Ritzerveld|2005 1. In practice, the effective absorption co- 
efficients will depend on the spectrum of the star and the 
recombination continuum, and their modification due to in- 
tervening absorption processes: the uncertainty is enhanced 
by the sharp dependence of the absorption cross section on 



the photon frequency above threshold (see Figure 2.2 in Os- 
terbrock and Ferland|2006 1. 

The calculations are normalized to a star emitting ion- 
izing photons at a rate 10'** s~^, generating an equilibrium 
Stromgren sphere of radius rs = 1 pc in all cases. This corre- 
sponds to a density of nu ~ ric = 180 cm~^ for the uniform 
density case. 

For direct comparison to Ritzerveld's work, we consider 
three main density distributions: uniform, and proportional 
to or beyond an empty core of radius 10~^ x rs and 
0.05 X rs, respectively. 



3 SOLUTION METHODS 

3.1 OTS and the diffusion limit 

The simplest approximation widely used in astrophysics is 
the case B approximation. Here, the diffuse field is assumed 
to be absorbed at the same point as it is emitted, so equa- 
tion ([2| is replaced by the on-the-spot (OTS) approximation 

n(l — x)ai/dif — '^n'x^ . (6) 
47r 

Note that the diffuse radiation intensity is isotropic in this 
approximation. With this assumption, equation Q simpli- 
fies to 



a-BU X 



and so 
dLdir/dr : 



i(l - a:)aoLdir 
47rr2 ' 



.2 2 2 

-47rr ctBTL X . 



(7) 



(8) 



Since the hydrogen ionization fraction is close to unity 
through most of an astrophysical Hll region, this equation 
may then be integrated to define the Stromgren radius, in 
agreement with equation ([5|. The ratio of the diffuse field 
photon density to that of the ionizing field is 



/dif,OTS = 



Ldir/47rr2 asai 



(9) 



which is roughly 60 per cent for equal absorption cross sec- 
tions, or ~ 10 per cent if we assume a lower absorption 
cross section for direct photons. The ratio of the diffuse 
fiux through a surface (of any orientation) to the direct fiux 
through a surface directed towards the star will be one quar- 
ter of this value (as in the analysis of |Cant6 et al.|1998[ who 
assume ao = ai, and thus find Fo = 0.15F*). 

If we wish to extend this analysis, we need to obtain a 
better estimate for Jdif = J IdifdO, to insert in equation (jsj). 
Writing equation ([2| as 

^ = + (10) 

where k = n(l — x)ai is the opacity and e — ain^x^ is the 
volume emissivity, then we can write 



Jdif = 



dn 



didii 
ds 



(11) 



This exact expression remains dependent on 7dif . If we apply 
the on-the-spot approximation Jdif — e/47rK to truncate the 
hierarchy at this approximation, the integration over angle is 
zero, which shows that the case B approximation is accurate 
to second order. 

To obtain more accurate results, however, we can sub- 
stitute for Jdif again using equation This gives 



Jdif 



d_ 

ds 



1 dJdi 



ds 



' dn. 

AtvkJ 



(12) 



Again, the second term in the integral will be zero. Assuming 
the remaining Jdif term is also isotropic to first order, it will 
be a function of radius alone, so 



d/dif 
ds 

d_ ( 1 d/dif 
ds 

and so 



/ 1 d/dif \ 
\ K ds J 



^ d/dif 
dr 

dr 



1 d/dif 



Hi dr 



+ sin 



J 1 d/dif 
Kr dr 



(13) 
(14) 
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J. ^ e ^ 1 d r dJdif 
K nr^ dr 3k dr ' 



(15) 



corresponding to Eddington approximation J = 3K (e.g. 
[Mihalas fc Weibel-Mihalas|1999| ). Thus given a distribution 
of ionization, we can solve for the diffuse field using this 
expression, in which the diffusive transport corrects the OTS 
intensity. Note, that since a; ~ 1 in the centre of the H ll 
region, the free path 1/k may be large. 



3.2 Ritzerveld's outward-only approximation 



Ritzerveld ( 2005 I solves an outward-only system 



dLdir i 2 d—l cLdir 

— ; — = -^O-An r 

dr 



(16) 



^ 2 d—l idif I ^ d-1 2 fT7\ 

= -4QAn r — h^r cun , (17) 

dr cLdir + -Ldit 

where for slab symmetry d = 1, ^ = 1, for cylindrical sym- 
metry d = 2, 5 = 27r and for spherical symmetry d = 3, 
^ = 47r. Note that these equations are given here for the 
total ionization intensity integrated over the surface at the 
specified radius. Ritzerveld states that they apply to the 
density of photons per unit volume, but this is clearly not 
the case, as if all the recombination coefficients are zero, then 
Ldir ~ Ldii ~ const. This system approximates the trans- 
port in the limit of nearly complete ionization within the 
Hll region. The factor c = ao/ai allows for the preferential 
absorption of the diffuse photons. 



The results presented by Ritzerveld ( 2005 1 can be gen- 



eralized to the case of arbitrary ratios of absorption cross 



section. Adding equations (161 and (171, we find 



/r , T \ f 2 d-1 

-r (-t^dir + -Ldif) = — ^QB". r , 

dr 



and hence 

Ltot = idir + idif = Li, 



f 2 d-1 J 

t^ot^n r dr. 



(18) 



(19) 



A general (implicit) expression for the direct and dif- 
fuse fields can be derived in terms of the total field, inde- 
pendent of the density distribution or geometry. Dividing 



equation ( 18 1 by equation ( 16 1, we find 



dLt, 



dLdir 

and hence 



Ob 
cqa 



Ltot 

Ldir 



(20) 



QfAC — Qb 



aicL* 



— (1 — c)aBLdir 



.(21) 



Given Ltotir) from equation (191 above, equation (21 1 may 
then be solved for Ldir (explicitly in the case c = 1, to recover 
the analytic solutions which Ritzerveld presents), and then 



3.3 Full transfer 

We use an iterative technique to find solutions, using similar 



techniques to those applied elsewhere (Hummer & Seaton 
1963[|Rubin|1968[ ). Here, we integrate the diffuse field using 
a large number of coaxial zones. 

We first derive an approximate ionization structure us- 
ing the case B approximation, and use this to estimate the 



outward-going diffuse field components. We integrate the 
diffuse field transport inwards through this structure. As 
part of this sweep, we calculate an improved estimate for 
the ionization structure taking account of the new inward 
radiation field beams. We then integrate the direct radia- 
tion and outward-going diffuse field outward through the 
resulting structure, using the inward diffuse field at its clos- 
est point to the star as the initial condition for the outward 
beams, updating the ionization structure again. We iterate 
these inward and outward sweeps until the ionization struc- 
ture has converged - this convergence is rapid in practice. 

The results for the ionization structure are very sensi- 
tive to the details of the numerical scheme for the case where 
the density varies as r~^, in particular the Stromgren radius 



can vary significantly. As noted by Franco et al. ( 1990 1, the 



ionization front can in this case escape to infinity: having the 
ionization front at 20 times the inner radius of the density 
distribution requires a rather precisely tuned value of the 
incident radiation field. This feeds through to the numerical 
sensitivity. 

As a result of this, considerable care has to be taken in 
the radiation field integration to ensure the rate at which 
photons are removed from the radiation field is consistent 
with the ionization rate within the zo ne ([Abel et al" 



1999 



WiUiamsl [20021 [Mellema et~aL] [20061 jWhalen fc Norman] 



20061. To do this, we will adapt the approach described by 



Williams ( 2002 1 to the case of multiple incident beams of 



radiation and diffuse sources. 

We specify the radiation field in coaxial cylindrical 
zones, which intersect spherical shells within which the phys- 
ical variables are assumed to be uniform. We take the i-th 
component of the density field to occupy a spherical shell of 
outer radius = lA, and the j-th group of diffuse photons 
to be those with an impact parameter between Rj-\ and 
Rj, where Rj = jA. 

Integrating equation ([sjl over a spatial zone, we have 

n(l - 



ffloLdir 

Anr^ 



+ 



ai 



Idit{r,n')dn' dV, 



where 

AT/ 47r / 3 3 



(22) 



(23) 



Substituting from the transport equations, ([TJ and ([2|, then 



asn^ AVi 



'^^I^ dQ'Anr^ dr, 
ds 



(24) 



where the problem is uniform in space, so the angular inte- 
gral over the zone volume is trivial. Combining the spatial 
integral r^dr with the beam solid angle integral dQ,' , how- 
ever, yields an integral which may be taken as over the zone 
volume for beams which are incident from a single direc- 
tion. This may be discretized as the difference in diffuse fiux 
within each radiation bin across the zone, 



d/di 



ds 



- df2'47rr^ dr = 



d/di 



ds 



-{b,s)ATT-2TTbdbds, (25) 
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where the integral may now be approximated as the change 
in the diffuse flux within the cylinder 

ImijAAj = / hii2Tvbdb (26) 

d/dif 



1/ 47r / 2 



R 



2\3/2 



The area of the cylindrical shell is 



(32) 



(33) 



ds 



-ds — /dif,i — /dif.s- 



(27) 



as the photons pass across the radial zone. Again, this has 
the interpretation that the difference between ionizing pho- 
tons into and out of the zone is balanced with the rate of 
recombinations in the zone which destroy such photons. 

In order to evaluate the differences in the direct and 
diffuse fluxes on either side of the zone, we discretize the 
transport equations. There are a range of choices for this 
discretization, which give different forms for the equations 



Hence the average path length of photons within the cylin- 
drical shell j through one of the the intersections with the 
spherical shell is 



(34) 



to be solved to evaluate the ionization fraction (Williams 



2002). 



As the source function from the recombinations means 
that there is no simple closed form solution to the ionization 
equilibrium problem in the zone, we instead choose a scheme 
where the outgoing beams are guaranteed to be positive. We 
integrate the transfer equation ([2| across a uniform zone to 
give 

/(s + As) ~ /(s)exp[-n(l - a::)aiAs] + 

Siit{x) (1 - exp [-naiil - x)As]) , (28) 

where 

„2 



5'dif(a;) 



ainx 



47rai(l — x) 



(29) 



is the diffuse field source function. The change in ionizing 
flux across the zone required for the ionization balance equa- 



tion in the form ( 24 1 is then 



7(s -t- As) - I{s) 



[Sdif(x) -J(s)] X 

(l-exp[-n(l-2:)aiAs]). (30) 



This expression must then be integrated over the cylindrical 
zone of impact parameters and the solid angle of photon 
directions to obtain a finite difference in the form needed to 



the first-order assumption that a single mean path length 
applies to all photons within the photon zone as they cross 
the radial zone. 

In order to perform the discrete integrals, we need to 
determine the volume of intersection between the radiation 
and spatial zones. The volume between these cylindrical and 
spherical shells, each of thickness A, is contained in two dis- 
joint components, apart from the case i — j. The radiation 
field in the latter case can still be advanced by the same 
formalism by splitting the zone in two at the intersection 
of the radiation beam with the equatorial plane. Note, how- 
ever, that the additional radiation flux values introduced 
into the diffuse fleld arrays as a result must not be included 
when calculating radiation field moments at the radial zone 
interfaces. 

Each component has volume 

Vi,3 = ^ {Vi-l,j + V^,j-l - Vi,j - Vi-l,j-l) , (31) 

where Vi,j is the volume within a coaxial sphere of radius Vi 
but outside a cylinder of radius Rj < ri, 



We use these distances as the chord lengths for propagation 
of the diffuse radiation propagation in the cylindrical shell 
labelled j through its intersection with the spherical shell 
labelled i in the ionization equilibrium solution, in order to 
obtain volume consistency. The resulting expression for the 
zone ionization, including all direct and diffuse beams, is 
solved using Newton-Raphson iteration. 

In all cases, we verify that at convergence the number 
of case B recombinations within the modelled H ll region 
balances the incoming ionizing flux to high accuracy. 



4 RESULTS 

In this section, we present the results of the diffuse radiation 
solutions using a variety of numerical schemes, compared to 
the full transfer solutions. The full transfer solutions have 
been discretized on a grid with spatial resolution of A = 
O.OOlrs out to 2rs. 

In Figure [l] we show the ionization fraction structures 
for the OTS and exact transfer solutions. These are seen to 
be very close to each other in all cases, excpet for close to the 
ionization front, particularly in the case of equal absorption 
coefficients. In the case of an r^^ density distribution, the 
outer edge of the H ll region is significantly broader. This is 
true for both the OTS and the full transfer solutions. The 
structure of the outer region of the ionized bubble in the 
case is highly sensitive to details of the internal structure, 
as can be seen from the relation 



be used in equation (241. On the right hand side, we make L* = AnoBn'^r'^irc ^ 



(35) 



so discretization errors at the 1 per cent level lead to changes 
in the outer radius of the ionized bubble rs /rc times greater, 
i.e. 20 per cent for the parameters chosen by IRitzerveld] 
(20051. In reality, it is likely that the ionized region would 



either be trapped in the core of the density distribution or 
extend to infinity; if the static solution were to be found in 
a region of such rapidly decreasing density, it is likely that 
the high sensitivity in its position would result in dynamical 
instability. 

In Figure [2] we show the direct and diffuse radiation in- 
tensities. As expected, the diffuse field intensity is higher in 
the case of equal absorption coefficients. Overall, the com- 
plete integration gives results roughly intermediate between 
the OTS and Ritzerveld approximations. 

The OTS approximation overestimates the diffuse field 
intensity in the core of the nebula, while the Ritzerveld 
approximation underestimates this. The direct field in the 
outer part of the nebula is reduced in the outer parts of the 
nebula in the Ritzerveld approximation, primarily as a re- 
sult of the lower ionization in the core. For the case of r~^. 
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Figure 1. Hll region structures. Left column, witli equal absorption coefficients, right column with diffuse coefficient 6 times direct. Top 
row uniform density, second row outside a small core, last row outside 0.05 core. Solid line is ionization fraction for full transfer 
solution, dotted is ionization fraction for OTS solution — in almost all cases, the ionization fractions are essentially identical. The lower 
panel in the plots shows Xf^n — xqtSi with the scale set to show differences internal to the Hll region. 
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Figure 2. Hll region structures. Left column, witli equal absorption coefficients, right column with diffuse coefficient 6 times direct. Top 
row uniform density, second row outside a small core, last row outside 0.05 core. Solid is direct photon rate, dashed is Anr^ 
times the local diffuse intensity. OTS results have square markers, the approximation of |Ritzerveld| p005[ l has triangular markers; the 
detailed transfer solution is plain. In general, the exact results lie between the two approximations; in lower right plot, the OTS and 
exact results are closely coincident. 
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Figure 3. Polar diagrams for the radiation at 0.1 (solid), 0.5 (dashed), 0.8 (dot-dashed) and 0.9 (dotted) xrg. Left column, with equal 
absorption coefficients, right column with diffuse absorption coefficient 6 times that for direct photons. Top row uniform density, second 
row outside a small core, last row outside 0.05 core. 
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Figure 4. Diffuse field components as a fraction of tlie direct field for the full-transfer solution. Left column, with equal absorption 
coefficients, right column with diffuse coefficient 6 times direct. Top row uniform density, second row outside a small core, last row 
outside 0.05 core. Solid is the total diffuse energy density, dashed is the flux incident tangentially, dot-dashed is the outward fiux, 
dotted is the inward fiux. Square and triangular markers on solid lines show the total diffuse energy density estimated by the OTS and 
Ritzerveld approximations, respectively. 
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Figure 5. Hll region structure and polar diagram for a thin shell with zero density out to 90 per cent of the Stromgren radius and 
uniform density beyond. Left column, with equal absorption coefficients, right column with diffuse absorption coefScient 6 times direct. 
The line style and colours of the curves have the same meanings as defined in previous figures. 
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the assumption that x = 1 within the Stromgren sphere 
means that this approximation also does not capture the 
broadening of the ionization front. 

In Figure |3j we show the angular distribution of the dif- 
fuse radiation at various positions within the nebula. Where 
the absorption coefficients of direct and diffuse radiation are 
the same, the diffuse radiation field is also beamed strongly 
forward through much of the nebula, because it is princi- 
pally emitted in the dense gas in the core. Hence, although 
the radiation energy may be distributed with the form char- 
acteristic of recombination emission, effects of a directed ra- 
diation field such as shadowing should still be in evidence. 
For the uniform density case with diffuse absorption coeffi- 
cient six time larger than direct, the radiation field is most 
strongly beamed at intermediate radii within the nebula: the 
directionality is small at r = 0.1 pc for both uniform density 
cases as this is close to the centre of symmetry, where as for 
the higher absorption case the diffuse field is also symmetric 
near to the edge of the nebula because the rdif ~ 1 surface 
samples a relatively small degree of anisotropy in the direct 
field. 

For the innermost radius of the distribution with 
the higher absorption coefficient for the direct field, the inner 
edge of the photoionized gas can be clearly seen in the diffuse 
radiation field. 

In Figure |4] we plot the diffuse radiation intensity as 
a fraction of the direct radiation intensity as a function 
of radius for our standard set of parameters. These results 
complement the polar diagrams at specific radii shown in 
Figure |3] As we have seen, in the case of equal absorption 
coefficients, the diffuse intensity can be higher than the di- 
rect intensity. However the ratio of the lateral diffuse flux 
to the direct flux, relevant to the illumination of the tails of 
cometary globules, reaches at most 60 per cent at the edge 
of a density distribution, and through most of the neb- 
ula is closer to the 15 per cent figure derived by |Cant6 et] 
al. (19981. In the case where the direct flux is less strongly 



absorbed, all the diffuse flux components are close to 2.5 
per cent of the direct flux through almost all of the region; 
the outward-going flux shows the most structure, but only 
reaches at most around 5.5 per cent of the direct flux at its 
maximum. 

This plot also shows the ratio of the diffuse to direct 
field in the OTS and Ritzerveld approximations. It is clear 
that the OTS approximation predicts a reasonable average 
value for this ratio, but overestimates the relative impor- 
tance of the diffuse field in the core of the nebula and 
underestimates it at the edge of the nebula. In contrast, 
Ritzerveld's formalism underestimates relative importance 
of the diffuse field in the core but overestimates it at the edge 
of the nebula, particularly for the case of equal absorption 
coefficients and steeply decreasing density distibutions (as 
we saw in Figure [2j this is primarily due to underestimating 
the direct field here). In both cases, the first approximation 
would be that the lateral diffuse field was one quarter the ra- 
diation intensity, which would be a significant overestimate 
for the case of equal absorption coefficients. 

In Figure |5] we show the same results as previously for 
a thin shell of uniform density between 0.9rs and rs- In this 
case, Ritzerveld's approximation obviously underestimates 
the diffuse field dramatically in the empty core of the neb- 
ula, as there are no recombinations in this region, while the 



OTS approximation overestimates the diffuse fiux. The ra- 
diation field is symmetric in angle at all our sample points, 
which are at radii less than the inner edge of the shell: this 
is to be expected, as the angle to the surface from oppo- 
sitely directed beams is the same and the surface brightness 
is constant. Particularly in the case with higher direct ab- 
sorption coefficient, limb-brightening is seen for the sample 
point at the inner surface of the shell, corresponding to the 
higher direct fiux (and hence source function) at the rdif ~ 1 
surface in this direction. 



4.1 Transfer using the diffusion approximation 

In Figure |6j we compare the accuracy of a diffusion approx- 
imation solution to the diffuse field intensity to that of the 
on-the-spot approximation. The diffusion results are almost 
everywhere more accurate than the OTS approximation, be- 
ing in error by less than 30 per cent even in the most diffi- 
cult, central region of the flow. Both schemes are accurate 
at the edge of the nebula, where the most interesting dy- 
namical structures occur, although the diffusion results are 
still significantly more precise. 

A diffusion approximation solver could be embedded in 
a dynamical code by updating the ionization solution con- 
sistently with the sources derived from the diffusion solver 
and iterating, in a similar fashion to that which we use for 
the full field transfer. Experimenting with updating the ion- 
ization structure after both inward and outward transfer are 
complete suggests that the solution converges more slowly 
than if the ionization structure is updated on each sweep. 
However, in the context of a dynamical scheme, the initial 
approximation for the ionization structure will have come 
from the previous timestep, and hence would be expected 
to not be far from the correct solution. 

The diffusion approximation allows processes such as 
the lateral illumination of the tails of photoionized globules 



to be modelled (Canto et al. 1998 Pavlakis et al. 20011 



However, the lateral illumination of shadowed tails may be 
overstimated in cases where the beaming of the diffuse field 



is expected to be strong (Mellema et al. 20061: Eddington 



tensor approaches (e.g. Gonzalez & Audit 2005 1 may im- 



prove the results in such limits. 



5 DISCUSSION 

5.1 Accuracy of the OTS approximation 

The accuracy of the OTS approximation may be evaluated 
by comparing the scale of variation of the diffuse field source 
function, H, to the mean free path for diffuse photons, A — 
l/n(l — x)ai. 

Using equation ( |29[ ) for the source function, and the 
OTS approximation equation ([7| for the ionization balance 

1 Qiao Ldir 

(36) 



4% aBfli 4:Tvr^ 



Note that we can still use the OTS approximation to close 
our equations at this higher order. The source function varies 
radially, so that 

1 dSdif 



H 



Sdii dr 



(37) 
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Figure 6. Accuracy of diffuse field determination for uniform density field. Left column, with equal absorption coefficients, right column 
with diffuse absorption coefficient 6 times that for direct photons. Solid line - ionization fraction; dotted line - OTS diffuse field divided 
by detailed transport diffuse field; dashed line - diffusion approximation diffuse field divided by detailed transport diffuse field. 




r/r„ 



Figure 7. Functional dependence from equation ( |43| for power- 
law density indexes m = (solid), m = 1 (dashed), m = 2 (dot- 
dashed). This function shows where the OTS approximation for 
the diffuse field is expected to be accurate (roughly for radii where 
the ordinate is less than ai/ao). 



A 
H 



ao 
ai 3 



2m 



(2(r/rs)= 



1 - 2m 



(43) 



(so long as m 7^ 3/2). In Figure [t] we show the form of the 
r.h.s. of equation (43 1 for ao — ai. This suggests that OTS 



approximation wiU be most accurate close to the outer edge 
of Hll regions. When ao = ai, the OTS approximation for 
the difTuse field is only marginally valid at any point in the 
nebula. For ao < ai, the OTS approximation may be valid 
for a larger region at the outside of the nebula. This is in 
accord with the results we have already seen in Figure |4] 
equation ( 43 1 is a simple formula which allows the accuracy 



to be assessed for a variety of Hll region structures. 

In conclusion, while the OTS approximation gives ac- 
curate results for the ionization structure in the H ll region, 
as we have seen in Figure [TJ the form of the diffuse field 
opacity n(l — x)ai emphasizes differences in the interior of 
the region, where x ~ 1. Hence, as demonstrated for certain 
H II region structures in Figure |4] and confirmed in general 
by the analysis in this section, the OTS approximation will 
give poor results for the diffuse radiation field, particularly 
when ao ~ ai. 



1 didir 2 

r 



dir dr 
Ai^r^ a^n? 2 

r ^ ~ ' 



(38) 
(39) 



where we have used equation (|8|. Now through most of the 
Hll region, a; ~ 1, and hence 



Airr'^a^n'^x^ a, 2 2> 

Ldir — — r-, N ^ — 47rr a^n A, 

n(i — 3;)ao ao 

and so 

A ft : LidiT 

H ~ 



ai 
ao 
ai 



1 + 
1 + 



2 l-rs 



(40) 

(41) 
(42) 



where the integral is the result of equation (|8|. For a density 
distribution n oc r"™ outside some specific radius. 



5.2 Dust effects on the diffuse field 

As has been mentioned, there is some observational evidence 
for dust absorption in well-observed H 11 regions. The effect 
of such dust absorption on nebular structure and emission 



has been widely studied, recently by [Arthur et al. (20041; 
Dopita et a"L] ([2006). 

Dust absorption will act to further limit the importance 
of the diffuse field. This is true independent of the relative 
dust absorption cross section at different frequencies, since 
the absorption of direct photons reduces the net source of 
diffuse photons. 

To obtain a more accurate estimate of the effects of 
dust, we assume a dust cross-section in EUV of ~ 
1.2 X 10~^^ cm^ per nucleon, with albedo lo ~ 0.5 (Bertoldi 
|fc Draine||1996[ ). The dust scattering function is primarily 
forward-beamed, so in the present analysis we will ignore 
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the dust-scattered continuum and take the effective absorp- 
tion cross-section to be a'^ — {1 ~ cj)(Tg. 

If we consider the OTS approximation, then for the dif- 
fuse field 



Jdif = 4:nS - 



ainx 



(1 — x)ai + (jg 



(44) 



We also have photoionization equilibrium, for which equa- 
tion ^ still holds. Substituting the OTS approximation ( 44 1 
into equation ([3|, we find 



Jdi 



(45) 



aAai[l + b{x)] 

where b{x) = a'^/ {l—x)ai. Rearranging, using qb ~ UA — ai, 
gives 

Jdif _ flp Qi 
Jdir ai QB + b{x)aA ' 



(46) 



Since in the presence of dust b{x) > 0, the dust absorption 
will always act to reduce the diffuse field, and will be signif- 



icant effect only where 6 ~ 1, i.e. (1 — a;) < fg/ai 
For a uniform nebula in the case B approximation. 



10" 



(47) 



so, if there is a sufficient dust optical depth, the effects of 
dust will be most significant close to the core of an H ll 
region. 



6 CONCLUSIONS 

We present detailed calculations of radiation transfer within 
a photoionized nebula, for a simplified physical model. We 
find that the diffuse field is indeed enhanced in the outer 
parts of such a nebula if there are steep density gradients in 
the Hll region, particularly if the absorption cross section 
for diffuse and direct photons is comparable. These circum- 
stances are not likely to be typical. We also find that when 
the diffuse radiation is relatively strong, it is also strongly 
beamed in the radial direction, and so the dynamical effect 
of the radiation field will in any case be similar to the direct 
illumination. 

We discuss the possibility of using a diffusion approx- 
imation for the diffuse radiation field in two- and three- 
dimensional radiation hydrodynamic calculations of H ll re- 
gions. This may allow the diffuse field effects to be calcu- 
lated to reasonable accuracy without requiring a full radia- 
tion transfer calculation. 

This paper has not considered in detail a number of ad- 
ditional processes which may effect the radiation field within 
the nebulae, such as the spectral hardening of radiation close 
to the ionization front, dust absorption and scattering, and 
the effects of helium and other heavy elements. This is left 
for future work. 
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